%main = getCoordination(main)
%   This function adds the coordination per element to the field
%   'coordination' by taking in account the amount of bound hydrogen atoms,
%   other bound atoms, bondtypes, charge and the expected coordination for
%   the "organic subset" of B,C,N,O,P,S,F,Cl,Br and I.

%%

function main = getCoordination(main)
%% Initialisation

%Loading the cell 'baseCoor' that contains the expected coordination for
%the members of the "organic subset".
load baseCoor.mat

%Defining the two string that will be using for error reporting. In this
%case that has been done here to make it easier to edit them and keeping
%the code easier to read.
str1 = [
    'The coordination of element %d in "main" was found to be greater\n'...
    'than any expected coordination. The program proceeded, but this \n'...
    'should be looked at. The maximum expected coordination was %d \n'...
    'but a coordination of %d was found.\n'];
str2 = [
    'The coordination was corrected for element %d in "main" to an \n'...
    'expected coordination while hydrogen were inputted for this \n'...
    'element. The used coordination (%d) therefore does not \n'...
    'corrospond with the inputted coordination (%d).\n'];

%% Main

%For all elements in 'main'
for j = 1:length(main)
    %Since this function needs the field 'hydrogen' to contain zeros
    %instead of empty spaces where not filled in this ifstatement is
    %required:
    if isempty(main(j).hydrogen)
        main(j).hydrogen = 0;
    end
    
    %The expected coordination is retrieved from 'baseCoor'
    wanted = [baseCoor{find(ismember({baseCoor{:,1}},main(j).element)),2}];
    %'baseCoor' However contains only the data for the "organic subset", so
    %for any other atoms the array 'wanted' will be empty. The function
    %required a numeric value instead of an empty space, so it is filled if
    %this is to occur.
    if isempty(wanted)
        wanted = 0;
    end
    
    %The coordination is actually the sum of all bound atoms (or groups)
    %but the wanted coordination stored in baseCoor is actually a valence,
    %so the sum of the bonded hydrogens and bondtypes is taken instead.
    coor = main(j).hydrogen + sum(main(j).bondstype);
    
    %If there are hydrogens bonded to the current atom
    if main(j).hydrogen ~= 0
        %Store 'coor' in 'coorstore'
        coorstor = coor;
    end
    
    %If the current atom has no charge
    if isempty(main(j).charge)
        %If the coordination is higher than the maximum expected
        %coordination and there is a wanted coordination known
        if coor > max(wanted) && wanted ~= 0
            %Current coordination should not be different, but is strange!
            %This can occur in special cases and should be reported, but
            %not meddled with. It can be an input mistake or on purpose.
            beep
            fprintf(str1,j,max(wanted),coor)
        elseif coor < min(wanted)
            %If the calculated coordination is smaller than the expected
            %coordination (because implicit hydrogen aren't accounted for)
            %change the coordination to be the minimum value in 'wanted'.
            coor = min(wanted);
        elseif length(wanted) > 1
            %If wanted contains multiple values
            %For every value in 'wanted' but the last
            for k = 1:length(wanted)-1
                %If the calculated coordination lies in between two values
                %in 'wanted'
                if coor > wanted(k) && coor < wanted(k+1)
                    %Then the coordination is the upper of the two values.
                    coor = wanted(k+1);
                end
            end
        else
            %'coor' Equals one of the wanted coordinations and is left alone
            %because of it.
        end
        
        %If the corrected coordination is not equal to the one calculated
        %with the hydrogens that are bound to it, print the second string
        %as an "error".
        if exist('coorstor','var') == 1
            if coorstor ~= coor
                beep
                fprintf(str2,j,coor,coorstor)
                clear coorstor
            end
        end
    else
        %When a charge is on the element, the coordination does not follow
        %the rules of expected coordination, since that only counts for
        %uncharged parts.
    end
    %Compenstate for double,triple and quadruple bonds
    coor = coor + length(main(j).bondstype) - sum(main(j).bondstype);
    %Store the coordination in the correct position in 'main'.
    main(j).coordination = coor;
end

%% End
end